Here, we will save DESeq2 objects for comparisons within species as well as comparisons across species.
library(tidyverse)
## ── Attaching core tidyverse packages ──────────────────────── tidyverse 2.0.0 ──
## ✔ dplyr 1.1.3 ✔ readr 2.1.4
## ✔ forcats 1.0.0 ✔ stringr 1.5.0
## ✔ ggplot2 3.4.4 ✔ tibble 3.2.1
## ✔ lubridate 1.9.3 ✔ tidyr 1.3.0
## ✔ purrr 1.0.2
## ── Conflicts ────────────────────────────────────────── tidyverse_conflicts() ──
## ✖ dplyr::filter() masks stats::filter()
## ✖ dplyr::lag() masks stats::lag()
## ℹ Use the conflicted package (<http://conflicted.r-lib.org/>) to force all conflicts to become errors
library(myTAI)
library(DESeq2)
## Loading required package: S4Vectors
## Loading required package: stats4
## Loading required package: BiocGenerics
##
## Attaching package: 'BiocGenerics'
##
## The following objects are masked from 'package:lubridate':
##
## intersect, setdiff, union
##
## The following objects are masked from 'package:dplyr':
##
## combine, intersect, setdiff, union
##
## The following objects are masked from 'package:stats':
##
## IQR, mad, sd, var, xtabs
##
## The following objects are masked from 'package:base':
##
## anyDuplicated, aperm, append, as.data.frame, basename, cbind,
## colnames, dirname, do.call, duplicated, eval, evalq, Filter, Find,
## get, grep, grepl, intersect, is.unsorted, lapply, Map, mapply,
## match, mget, order, paste, pmax, pmax.int, pmin, pmin.int,
## Position, rank, rbind, Reduce, rownames, sapply, setdiff, sort,
## table, tapply, union, unique, unsplit, which.max, which.min
##
##
## Attaching package: 'S4Vectors'
##
## The following objects are masked from 'package:lubridate':
##
## second, second<-
##
## The following objects are masked from 'package:dplyr':
##
## first, rename
##
## The following object is masked from 'package:tidyr':
##
## expand
##
## The following objects are masked from 'package:base':
##
## expand.grid, I, unname
##
## Loading required package: IRanges
##
## Attaching package: 'IRanges'
##
## The following object is masked from 'package:lubridate':
##
## %within%
##
## The following objects are masked from 'package:dplyr':
##
## collapse, desc, slice
##
## The following object is masked from 'package:purrr':
##
## reduce
##
## Loading required package: GenomicRanges
## Loading required package: GenomeInfoDb
## Loading required package: SummarizedExperiment
## Loading required package: MatrixGenerics
## Loading required package: matrixStats
##
## Attaching package: 'matrixStats'
##
## The following object is masked from 'package:dplyr':
##
## count
##
##
## Attaching package: 'MatrixGenerics'
##
## The following objects are masked from 'package:matrixStats':
##
## colAlls, colAnyNAs, colAnys, colAvgsPerRowSet, colCollapse,
## colCounts, colCummaxs, colCummins, colCumprods, colCumsums,
## colDiffs, colIQRDiffs, colIQRs, colLogSumExps, colMadDiffs,
## colMads, colMaxs, colMeans2, colMedians, colMins, colOrderStats,
## colProds, colQuantiles, colRanges, colRanks, colSdDiffs, colSds,
## colSums2, colTabulates, colVarDiffs, colVars, colWeightedMads,
## colWeightedMeans, colWeightedMedians, colWeightedSds,
## colWeightedVars, rowAlls, rowAnyNAs, rowAnys, rowAvgsPerColSet,
## rowCollapse, rowCounts, rowCummaxs, rowCummins, rowCumprods,
## rowCumsums, rowDiffs, rowIQRDiffs, rowIQRs, rowLogSumExps,
## rowMadDiffs, rowMads, rowMaxs, rowMeans2, rowMedians, rowMins,
## rowOrderStats, rowProds, rowQuantiles, rowRanges, rowRanks,
## rowSdDiffs, rowSds, rowSums2, rowTabulates, rowVarDiffs, rowVars,
## rowWeightedMads, rowWeightedMeans, rowWeightedMedians,
## rowWeightedSds, rowWeightedVars
##
## Loading required package: Biobase
## Welcome to Bioconductor
##
## Vignettes contain introductory material; view with
## 'browseVignettes()'. To cite Bioconductor, see
## 'citation("Biobase")', and for packages 'citation("pkgname")'.
##
##
## Attaching package: 'Biobase'
##
## The following object is masked from 'package:MatrixGenerics':
##
## rowMedians
##
## The following objects are masked from 'package:matrixStats':
##
## anyMissing, rowMedians
library(IHW)
##
## Attaching package: 'IHW'
##
## The following object is masked from 'package:ggplot2':
##
## alpha
sample_info_ec <-
readr::read_csv(file = "data/sample_info_ec.csv")
## Rows: 57 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (9): SampleName, SeqProt, Species, stage_tissue, Stage, Tissue, Experime...
## dbl (1): Replicate
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
sample_info_fd <-
readr::read_csv(file = "data/sample_info_fd.csv")
## Rows: 34 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (9): SampleName, SeqProt, Species, stage_tissue, Stage, Tissue, Experime...
## dbl (1): Replicate
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
sample_info_fs <-
readr::read_csv(file = "data/sample_info_fs.csv")
## Rows: 48 Columns: 10
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (9): SampleName, SeqProt, Species, stage_tissue, Stage, Tissue, Experime...
## dbl (1): Replicate
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_counts <-
readr::read_csv("data/Ec_counts.csv")
## Rows: 11582 Columns: 58
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (57): M_meiospore_low_1, M_meiospore_low_2, M_meiospore_low_3, F_meiospo...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_counts <-
readr::read_csv("data/Fd_counts.csv")
## Rows: 7907 Columns: 35
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (34): Fd_f_gametes_1, Fd_f_gametes_2, Fd_f_gametes_3, Fd_m_gametes_1, Fd...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_counts <-
readr::read_csv("data/Fs_counts.csv")
## Rows: 8291 Columns: 49
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (48): Fs_F_gametes_1, Fs_F_gametes_2, Fs_F_gametes_3, Fs_M_gametes_1, Fs...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Ec_counts.denoised <-
readr::read_csv("data/Ec_counts_denoised.csv")
## Rows: 2695 Columns: 58
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (57): M_meiospore_low_1, M_meiospore_low_2, M_meiospore_low_3, F_meiospo...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fd_counts.denoised <-
readr::read_csv("data/Fd_counts_denoised.csv")
## Rows: 2010 Columns: 35
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (34): Fd_f_gametes_1, Fd_f_gametes_2, Fd_f_gametes_3, Fd_m_gametes_1, Fd...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Fs_counts.denoised <-
readr::read_csv("data/Fs_counts_denoised.csv")
## Rows: 2214 Columns: 49
## ── Column specification ────────────────────────────────────────────────────────
## Delimiter: ","
## chr (1): GeneID
## dbl (48): Fs_F_gametes_1, Fs_F_gametes_2, Fs_F_gametes_3, Fs_M_gametes_1, Fs...
##
## ℹ Use `spec()` to retrieve the full column specification for this data.
## ℹ Specify the column types or set `show_col_types = FALSE` to quiet this message.
Function to reduce the number of lines of code.
# the count data is in tibble format as I have saved it.
create_dds_embyro <-
function(count_tibble, metatable, stage_names, rename_list){
# some error messages
if(!is_tibble(count_tibble)){return("is not tibble")}
if(!colnames(count_tibble)[1] == "GeneID"){return("make sure the first column is titled \"GeneID\"")}
# filter the metatable
metatable.filt <-
metatable %>%
dplyr::filter(Stage %in% stage_names) %>%
dplyr::arrange(base::match(Stage, stage_names))
# filter the counts to be loaded as a dds object
count.embryo <-
count_tibble %>%
dplyr::select(
"GeneID" |
metatable.filt$LibName # this would be different for Ectocarpus
) %>%
dplyr::rename_with(~rename_list[.x], -GeneID) %>%
tibble::remove_rownames() %>%
tibble::column_to_rownames(var="GeneID")
# creating a ddsObj for the embryo stages
simple.model <-
as.formula(~Stage)
ddsObj.raw <-
DESeqDataSetFromMatrix(
round(count.embryo,0),
colData = metatable.filt,
design = simple.model)
}
Fs_embryo_stage_names <-
c("24H", "48H", "1w", "3w", "4w")
sample_info_fs[13,2] <- 4
sample_info_fs[14,2] <- 5
sample_info_fs[15,2] <- 6
rename_list <-
sample_info_fs %>%
mutate(new_colname = paste0(Stage, "_", Sex, "_", Replicate)) %>%
select(LibName, new_colname) %>%
deframe()
Fs_stage_ddsObj.raw <-
create_dds_embyro(
count_tibble = Fs_counts,
metatable = sample_info_fs,
stage_names = Fs_embryo_stage_names,
rename_list = rename_list)
## converting counts to integer mode
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
rename_list <-
select(sample_info_fd, LibName, Stage, Sex, Replicate) %>%
mutate(new_colname = paste0(Stage, "_", Sex, "_", Replicate)) %>%
select(LibName, new_colname) %>%
deframe()
Fd_embryo_stage_names <-
c("E1", "E2", "E3", "E4", "E5", "E6")
Fd_stage_ddsObj.raw <-
create_dds_embyro(
count_tibble = Fd_counts,
metatable = sample_info_fd,
stage_names = Fd_embryo_stage_names,
rename_list = rename_list)
## converting counts to integer mode
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
# significantly DEG (LRT)
Fs_ddsObj <-
DESeq2::DESeq(Fs_stage_ddsObj.raw, test="LRT", reduced=~1)
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
Fs_resLRT <-
DESeq2::results(Fs_ddsObj)
Fd_ddsObj <-
DESeq2::DESeq(Fd_stage_ddsObj.raw, test="LRT", reduced=~1)
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
Fd_resLRT <-
DESeq2::results(Fd_ddsObj)
# saving with padj filters. With LRT, log foldchange thresholds do not matter.
# see https://hbctraining.github.io/DGE_workshop/lessons/08_DGE_LRT.html
padj_filt = 0.001
tibble::as_tibble(Fs_resLRT, rownames = "GeneID") %>%
dplyr::filter(padj < padj_filt)
tibble::as_tibble(Fd_resLRT, rownames = "GeneID") %>%
dplyr::filter(padj < padj_filt)
tibble::as_tibble(Fs_resLRT, rownames = "GeneID") %>%
ggplot2::ggplot(aes(padj)) +
ggplot2::geom_histogram() +
ggplot2::scale_y_log10() +
ggplot2::labs(
title = "LRT",
subtitle = "Fucus serratus",
x = "padj",
y = "density (log-scaled)"
)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 60 rows containing non-finite values (`stat_bin()`).
tibble::as_tibble(Fd_resLRT, rownames = "GeneID") %>%
ggplot2::ggplot(aes(padj)) +
ggplot2::geom_histogram() +
ggplot2::scale_y_log10() +
ggplot2::labs(
title = "LRT",
subtitle = "Fucus distichus",
x = "padj",
y = "density (log-scaled)"
)
## `stat_bin()` using `bins = 30`. Pick better value with `binwidth`.
## Warning: Removed 24 rows containing non-finite values (`stat_bin()`).
plot(
sort(-log10(Fs_resLRT$padj),decreasing = TRUE),
xlab="Genes",
ylab="-log10(p-value)",
main = "Why select 1000 genes?",
sub = "Fucus serratus")
abline(v=1000, col="blue")
plot(
sort(-log10(Fd_resLRT$padj),decreasing = TRUE),
xlab="Genes",
ylab="-log10(p-value)",
main = "Why select 1000 genes?",
sub = "Fucus distichus")
abline(v=1000, col="blue")
# select top 1000 DEGs instead.
Fs_LRT.1000 <-
tibble::as_tibble(Fs_resLRT, rownames = "GeneID") %>%
arrange(padj) %>%
head(n=1000)
# Max padj:3.828e-40
Fd_LRT.1000 <-
tibble::as_tibble(Fd_resLRT, rownames = "GeneID") %>%
arrange(padj) %>%
head(n=1000)
# Max padj:4.405e-37
Fs_LRT.1000 %>%
readr::write_csv("data/Fs_LRT_1000.csv")
Fd_LRT.1000 %>%
readr::write_csv("data/Fd_LRT_1000.csv")
6,271 F. serratus and 5,845 F. distichus genes pass
the padj filter. This is too much. I have filtered the top 1000 most
differentially expressed genes instead by ranking the padj according to
the LRT. With the LRT, log2FoldChange is
irrelevant.
Furthermore, I have to look a bit more into the LRT.
This is because I had thought the LRT tests for genes that
can be explained by modelling Stage vs NULL. However, the
log2 fold change (MLE): Stage E5 vs E1 makes me worried.
However, I think it is fine when ranking with p-values or
padj because that is calculated via
LRT p-value: '~ Stage' vs '~ 1', which is different to just
applying DESeq2::results() on Wald,
i.e. Wald test p-value: Stage 4w vs 1w.
Whew.
> Fs_resLRT
log2 fold change (MLE): Stage 4w vs 1w
LRT p-value: '~ Stage' vs '~ 1'
DataFrame with 8291 rows and 6 columns
> DESeq2::results(Fs_ddsObj.wald)
log2 fold change (MLE): Stage 4w vs 1w
Wald test p-value: Stage 4w vs 1w
DataFrame with 8291 rows and 6 columns
However, I think is fine.
This uses the Wald’s test.
summarise_contrasts <-
function(
ddsObj,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names,
IHW = FALSE,
output_nsig = TRUE){
contrast.res <-
list()
if(!IHW){
for (i in 1:(length(stage_names)-1)) {
contrast.res[[i]] <-
DESeq2::results(
ddsObj,
lfcThreshold = lfcThreshold,
altHypothesis = altHypothesis,
contrast = c("Stage", stage_names[i+1], stage_names[i]))
}
}
else if(IHW){
for (i in 1:(length(stage_names)-1)) {
contrast.res[[i]] <-
DESeq2::results(
ddsObj,
lfcThreshold = lfcThreshold,
altHypothesis = altHypothesis,
filterFun = ihw,
contrast = c("Stage", stage_names[i+1], stage_names[i]))
}
}
if(!output_nsig){
return(contrast.res)
}
contrast.sig.count <-
tibble::tibble(
"contrast" = deseq2_contrasts,
"n_sig" = 0,
"lfcThreshold" = lfcThreshold,
"padj_threshold" = padj_threshold,
"altHypothesis" = altHypothesis,
"IHW" = IHW
)
for (i in 1:length(contrast.res)) {
contrast.sig.count$n_sig[i] <-
count(contrast.res[[i]]$padj < padj_threshold, na.rm = T)
}
return(contrast.sig.count)
}
Fs_ddsObj.wald <-
DESeq2::DESeq(Fs_stage_ddsObj.raw, test="Wald")
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
Fd_ddsObj.wald <-
DESeq2::DESeq(Fd_stage_ddsObj.raw, test="Wald")
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## final dispersion estimates
## fitting model and testing
Fucus serratus
Fs_embryo_stage_names <-
c("24H", "48H", "1w", "3w", "4w")
deseq2_contrasts <- c("24H->48H", "48H->1w", "1w->3w", "3w->4w")
Fs_ddsObj.wald_LFC1 <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC2 <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC3 <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC1_IHW <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
IHW = TRUE)
Fs_ddsObj.wald_LFC2_IHW <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
IHW = TRUE)
Fs_ddsObj.wald_LFC3_IHW <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
IHW = TRUE)
Fs_ddsObj.wald_padj0.1 <-
base::rbind(Fs_ddsObj.wald_LFC1,
Fs_ddsObj.wald_LFC2,
Fs_ddsObj.wald_LFC3,
Fs_ddsObj.wald_LFC1_IHW,
Fs_ddsObj.wald_LFC2_IHW,
Fs_ddsObj.wald_LFC3_IHW)
Fs_ddsObj.wald_LFC1 <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC2 <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC3 <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC1_IHW <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
IHW = TRUE)
Fs_ddsObj.wald_LFC2_IHW <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
IHW = TRUE)
Fs_ddsObj.wald_LFC3_IHW <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
IHW = TRUE)
Fs_ddsObj.wald_padj0.05 <-
base::rbind(Fs_ddsObj.wald_LFC1,
Fs_ddsObj.wald_LFC2,
Fs_ddsObj.wald_LFC3,
Fs_ddsObj.wald_LFC1_IHW,
Fs_ddsObj.wald_LFC2_IHW,
Fs_ddsObj.wald_LFC3_IHW)
# significantly nonDEG
Fs_ddsObj.wald_LFC1_less <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC2_less <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC3_less <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_less_padj0.1 <-
base::rbind(Fs_ddsObj.wald_LFC1_less,
Fs_ddsObj.wald_LFC2_less,
Fs_ddsObj.wald_LFC3_less)
Fs_ddsObj.wald_LFC1_less <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC2_less <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_LFC3_less <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.wald_less_padj0.05 <-
base::rbind(Fs_ddsObj.wald_LFC1_less,
Fs_ddsObj.wald_LFC2_less,
Fs_ddsObj.wald_LFC3_less)
Fucus distichus
Fd_embryo_stage_names <-
c("E1", "E2", "E3", "E4", "E5", "E6")
deseq2_contrasts <- c("E1->E2", "E2->E3", "E3->E4", "E4->E5", "E5->E6")
Fd_ddsObj.wald_LFC1 <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC2 <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC3 <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC1_IHW <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
IHW = TRUE)
Fd_ddsObj.wald_LFC2_IHW <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
IHW = TRUE)
Fd_ddsObj.wald_LFC3_IHW <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
IHW = TRUE)
Fd_ddsObj.wald_padj0.1 <-
base::rbind(Fd_ddsObj.wald_LFC1,
Fd_ddsObj.wald_LFC2,
Fd_ddsObj.wald_LFC3,
Fd_ddsObj.wald_LFC1_IHW,
Fd_ddsObj.wald_LFC2_IHW,
Fd_ddsObj.wald_LFC3_IHW)
Fd_ddsObj.wald_LFC1 <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC2 <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC3 <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC1_IHW <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
IHW = TRUE)
Fd_ddsObj.wald_LFC2_IHW <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
IHW = TRUE)
Fd_ddsObj.wald_LFC3_IHW <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
IHW = TRUE)
Fd_ddsObj.wald_padj0.05 <-
base::rbind(Fd_ddsObj.wald_LFC1,
Fd_ddsObj.wald_LFC2,
Fd_ddsObj.wald_LFC3,
Fd_ddsObj.wald_LFC1_IHW,
Fd_ddsObj.wald_LFC2_IHW,
Fd_ddsObj.wald_LFC3_IHW)
# significantly nonDEG
Fd_ddsObj.wald_LFC1_less <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC2_less <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC3_less <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_less_padj0.1 <-
base::rbind(Fd_ddsObj.wald_LFC1_less,
Fd_ddsObj.wald_LFC2_less,
Fd_ddsObj.wald_LFC3_less)
Fd_ddsObj.wald_LFC1_less <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC2_less <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_LFC3_less <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.wald_less_padj0.05 <-
base::rbind(Fd_ddsObj.wald_LFC1_less,
Fd_ddsObj.wald_LFC2_less,
Fd_ddsObj.wald_LFC3_less)
# combine everything together
Fs_ddsObj.wald_all <-
base::rbind(Fs_ddsObj.wald_padj0.1,
Fs_ddsObj.wald_padj0.05,
Fs_ddsObj.wald_less_padj0.1,
Fs_ddsObj.wald_less_padj0.05)
Fs_ddsObj.wald_all$contrast <-
base::factor(Fs_ddsObj.wald_all$contrast, levels = unique(Fs_ddsObj.wald_all$contrast))
Fd_ddsObj.wald_all <-
base::rbind(Fd_ddsObj.wald_padj0.1,
Fd_ddsObj.wald_padj0.05,
Fd_ddsObj.wald_less_padj0.1,
Fd_ddsObj.wald_less_padj0.05)
Fd_ddsObj.wald_all$contrast <-
base::factor(Fd_ddsObj.wald_all$contrast, levels = unique(Fd_ddsObj.wald_all$contrast))
The difference between performing IHW or not is almost the same.
Fs_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "greaterAbs" & IHW == FALSE) %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10() +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus serratus",
subtitle = "DEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
## Warning: Transformation introduced infinite values in continuous y-axis
## Transformation introduced infinite values in continuous y-axis
Fs_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "greaterAbs" & IHW == TRUE) %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10() +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus serratus",
subtitle = "DEG (padj = ihw)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
## Warning: Transformation introduced infinite values in continuous y-axis
## Transformation introduced infinite values in continuous y-axis
Fs_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10(limits = c(1, 10000)) +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus serratus",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
Fs_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::ylim(0, 10000) +
ggplot2::labs(
y = "Number of genes",
x = "Developmental transition",
title = "Fucus serratus",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
Fd_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "greaterAbs" & IHW == FALSE) %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10() +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus distichus",
subtitle = "DEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
## Warning: Transformation introduced infinite values in continuous y-axis
## Transformation introduced infinite values in continuous y-axis
Fd_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "greaterAbs" & IHW == TRUE) %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10() +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus distichus",
subtitle = "DEG (padj = ihw)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
## Warning: Transformation introduced infinite values in continuous y-axis
## Transformation introduced infinite values in continuous y-axis
Fd_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10(limits = c(1, 10000)) +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus distichus",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
Fd_ddsObj.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::ylim(0, 10000) +
ggplot2::labs(
y = "Number of genes",
x = "Developmental transition",
title = "Fucus distichus",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
# Fs_embryo_stage_names <-
# c("24H", "48H", "1w", "3w", "4w")
# sample_info_fs[13,2] <- 4
# sample_info_fs[14,2] <- 5
# sample_info_fs[15,2] <- 6
rename_list <-
sample_info_fs %>%
mutate(new_colname = paste0(Stage, "_", Sex, "_", Replicate)) %>%
select(LibName, new_colname) %>%
deframe()
Fs_stage_ddsObj.denoised.raw <-
create_dds_embyro(
count_tibble = Fs_counts.denoised,
metatable = sample_info_fs,
stage_names = Fs_embryo_stage_names,
rename_list = rename_list)
## converting counts to integer mode
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
rename_list <-
select(sample_info_fd, LibName, Stage, Sex, Replicate) %>%
mutate(new_colname = paste0(Stage, "_", Sex, "_", Replicate)) %>%
select(LibName, new_colname) %>%
deframe()
# Fd_embryo_stage_names <-
# c("E1", "E2", "E3", "E4", "E5", "E6")
Fd_stage_ddsObj.denoised.raw <-
create_dds_embyro(
count_tibble = Fd_counts.denoised,
metatable = sample_info_fd,
stage_names = Fd_embryo_stage_names,
rename_list = rename_list)
## converting counts to integer mode
## Warning in DESeqDataSet(se, design = design, ignoreRank): some variables in
## design formula are characters, converting to factors
Fs_ddsObj.denoised.wald <-
DESeq2::DESeq(Fs_stage_ddsObj.denoised.raw, test="Wald")
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## -- note: fitType='parametric', but the dispersion trend was not well captured by the
## function: y = a/x + b, and a local regression fit was automatically substituted.
## specify fitType='local' or 'mean' to avoid this message next time.
## final dispersion estimates
## fitting model and testing
Fd_ddsObj.denoised.wald <-
DESeq2::DESeq(Fd_stage_ddsObj.denoised.raw, test="Wald")
## estimating size factors
## estimating dispersions
## gene-wise dispersion estimates
## mean-dispersion relationship
## -- note: fitType='parametric', but the dispersion trend was not well captured by the
## function: y = a/x + b, and a local regression fit was automatically substituted.
## specify fitType='local' or 'mean' to avoid this message next time.
## final dispersion estimates
## fitting model and testing
There is no need to use IHW because of the warning
Only 1 bin; IHW reduces to Benjamini Hochberg (uniform weights).
Fucus serratus
Fs_embryo_stage_names <-
c("24H", "48H", "1w", "3w", "4w")
deseq2_contrasts <- c("24H->48H", "48H->1w", "1w->3w", "3w->4w")
Fs_ddsObj.denoised.wald_LFC1 <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC2 <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC3 <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_padj0.1 <-
base::rbind(Fs_ddsObj.denoised.wald_LFC1,
Fs_ddsObj.denoised.wald_LFC2,
Fs_ddsObj.denoised.wald_LFC3)
Fs_ddsObj.denoised.wald_LFC1 <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC2 <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC3 <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_padj0.05 <-
base::rbind(Fs_ddsObj.denoised.wald_LFC1,
Fs_ddsObj.denoised.wald_LFC2,
Fs_ddsObj.denoised.wald_LFC3)
# significantly nonDEG
Fs_ddsObj.denoised.wald_LFC1_less <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC2_less <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC3_less <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_less_padj0.1 <-
base::rbind(Fs_ddsObj.denoised.wald_LFC1_less,
Fs_ddsObj.denoised.wald_LFC2_less,
Fs_ddsObj.denoised.wald_LFC3_less)
Fs_ddsObj.denoised.wald_LFC1_less <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC2_less <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_LFC3_less <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names)
Fs_ddsObj.denoised.wald_less_padj0.05 <-
base::rbind(Fs_ddsObj.denoised.wald_LFC1_less,
Fs_ddsObj.denoised.wald_LFC2_less,
Fs_ddsObj.denoised.wald_LFC3_less)
Fucus distichus
Fd_embryo_stage_names <-
c("E1", "E2", "E3", "E4", "E5", "E6")
deseq2_contrasts <- c("E1->E2", "E2->E3", "E3->E4", "E4->E5", "E5->E6")
Fd_ddsObj.denoised.wald_LFC1 <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC2 <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC3 <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_padj0.1 <-
base::rbind(Fd_ddsObj.denoised.wald_LFC1,
Fd_ddsObj.denoised.wald_LFC2,
Fd_ddsObj.denoised.wald_LFC3)
Fd_ddsObj.denoised.wald_LFC1 <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC2 <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC3 <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_padj0.05 <-
base::rbind(Fd_ddsObj.denoised.wald_LFC1,
Fd_ddsObj.denoised.wald_LFC2,
Fd_ddsObj.denoised.wald_LFC3)
# significantly nonDEG
Fd_ddsObj.denoised.wald_LFC1_less <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC2_less <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC3_less <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.1,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_less_padj0.1 <-
base::rbind(Fd_ddsObj.denoised.wald_LFC1_less,
Fd_ddsObj.denoised.wald_LFC2_less,
Fd_ddsObj.denoised.wald_LFC3_less)
Fd_ddsObj.denoised.wald_LFC1_less <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC2_less <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 2,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_LFC3_less <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 3,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names)
Fd_ddsObj.denoised.wald_less_padj0.05 <-
base::rbind(Fd_ddsObj.denoised.wald_LFC1_less,
Fd_ddsObj.denoised.wald_LFC2_less,
Fd_ddsObj.denoised.wald_LFC3_less)
# combine everything together
Fs_ddsObj.denoised.wald_all <-
base::rbind(Fs_ddsObj.denoised.wald_padj0.1,
Fs_ddsObj.denoised.wald_padj0.05,
Fs_ddsObj.denoised.wald_less_padj0.1,
Fs_ddsObj.denoised.wald_less_padj0.05)
Fs_ddsObj.denoised.wald_all$contrast <-
base::factor(Fs_ddsObj.denoised.wald_all$contrast, levels = unique(Fs_ddsObj.denoised.wald_all$contrast))
Fd_ddsObj.denoised.wald_all <-
base::rbind(Fd_ddsObj.denoised.wald_padj0.1,
Fd_ddsObj.denoised.wald_padj0.05,
Fd_ddsObj.denoised.wald_less_padj0.1,
Fd_ddsObj.denoised.wald_less_padj0.05)
Fd_ddsObj.denoised.wald_all$contrast <-
base::factor(Fd_ddsObj.denoised.wald_all$contrast, levels = unique(Fd_ddsObj.denoised.wald_all$contrast))
The difference between performing IHW or not is almost the same.
Fs_ddsObj.denoised.wald_all %>%
dplyr::filter(altHypothesis == "greaterAbs" & IHW == FALSE) %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10() +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus serratus [denoised]",
subtitle = "DEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
## Warning: Transformation introduced infinite values in continuous y-axis
## Transformation introduced infinite values in continuous y-axis
Fs_ddsObj.denoised.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10(limits = c(1, 10000)) +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus serratus [denoised]",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
Fs_ddsObj.denoised.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::ylim(0, 10000) +
ggplot2::labs(
y = "Number of genes",
x = "Developmental transition",
title = "Fucus serratus [denoised]",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
Fd_ddsObj.denoised.wald_all %>%
dplyr::filter(altHypothesis == "greaterAbs" & IHW == FALSE) %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10() +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus distichus [denoised]",
subtitle = "DEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
## Warning: Transformation introduced infinite values in continuous y-axis
## Transformation introduced infinite values in continuous y-axis
Fd_ddsObj.denoised.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_y_log10(limits = c(1, 10000)) +
# ggplot2::annotation_logticks(sides = "l") +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::labs(
y = "Number of genes (log-scaled)",
x = "Developmental transition",
title = "Fucus distichus [denoised]",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
Fd_ddsObj.denoised.wald_all %>%
dplyr::filter(altHypothesis == "lessAbs") %>%
ggplot2::ggplot(
aes(
x = contrast,
y = n_sig,
colour = as.factor(lfcThreshold),
group = as.factor(lfcThreshold)
)
) +
ggplot2::geom_line(linewidth = 2, alpha = 0.5) +
ggplot2::geom_point(size = 3, alpha = 0.5) +
ggplot2::facet_grid(rows = ggplot2::vars(padj_threshold), labeller = labeller(padj_threshold = label_both)) +
ggplot2::scale_x_discrete(guide = guide_axis(angle = 30)) +
ggplot2::ylim(0, 10000) +
ggplot2::labs(
y = "Number of genes",
x = "Developmental transition",
title = "Fucus distichus [denoised]",
subtitle = "nonDEG (padj = BH)"
) +
ggplot2::theme_bw() +
ggplot2::theme(legend.position="bottom") +
ggplot2::guides(colour=guide_legend(title="logFoldChange threshold")) +
ggsci::scale_color_npg()
The count isn’t very interesting. So using the same function I made
summarise_contrasts() with
output_nsig = FALSE. All the data is based on the following
thresholds lfcThreshold = 1 and
padj_threshold
deseq2_contrasts <- c("24H->48H", "48H->1w", "1w->3w", "3w->4w")
Fs_ddsObj.wald_DEGlist <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
output_nsig = FALSE)
Fs_ddsObj.wald_nonDEGlist <-
summarise_contrasts(
Fs_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names,
output_nsig = FALSE)
deseq2_contrasts <- c("E1->E2", "E2->E3", "E3->E4", "E4->E5", "E5->E6")
Fd_ddsObj.wald_DEGlist <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
output_nsig = FALSE)
Fd_ddsObj.wald_nonDEGlist <-
summarise_contrasts(
Fd_ddsObj.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names,
output_nsig = FALSE)
deseq2_contrasts <- c("24H->48H", "48H->1w", "1w->3w", "3w->4w")
Fs_ddsObj.denoised.wald_DEGlist <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fs_embryo_stage_names,
output_nsig = FALSE)
Fs_ddsObj.denoised.wald_nonDEGlist <-
summarise_contrasts(
Fs_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fs_embryo_stage_names,
output_nsig = FALSE)
deseq2_contrasts <- c("E1->E2", "E2->E3", "E3->E4", "E4->E5", "E5->E6")
Fd_ddsObj.denoised.wald_DEGlist <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "greaterAbs",
stage_names = Fd_embryo_stage_names,
output_nsig = FALSE)
Fd_ddsObj.denoised.wald_nonDEGlist <-
summarise_contrasts(
Fd_ddsObj.denoised.wald,
deseq2_contrasts,
lfcThreshold = 1,
padj_threshold = 0.05,
altHypothesis = "lessAbs",
stage_names = Fd_embryo_stage_names,
output_nsig = FALSE)
Fs_ddsObj.wald_DEGlist %>%
base::saveRDS("data/Fs_wald_DEG_list.rds")
Fd_ddsObj.wald_DEGlist %>%
base::saveRDS("data/Fd_wald_DEG_list.rds")
Fs_ddsObj.wald_nonDEGlist %>%
base::saveRDS("data/Fs_wald_nonDEG_list.rds")
Fd_ddsObj.wald_nonDEGlist %>%
base::saveRDS("data/Fd_wald_nonDEG_list.rds")
Fs_ddsObj.denoised.wald_DEGlist %>%
base::saveRDS("data/Fs_denoised_wald_DEG_list.rds")
Fd_ddsObj.denoised.wald_DEGlist %>%
base::saveRDS("data/Fd_denoised_wald_DEG_list.rds")
Fs_ddsObj.denoised.wald_nonDEGlist %>%
base::saveRDS("data/Fs_denoised_wald_nonDEG_list.rds")
Fd_ddsObj.denoised.wald_nonDEGlist %>%
base::saveRDS("data/Fd_denoised_wald_nonDEG_list.rds")
# # to load data
# Fs_ddsObj.wald_DEGlist <-
# base::readRDS("data/Fs_wald_DEG_list.rds")
# rbind(
# Fs_ddsObj.wald_DEGlist[[1]] %>%
# tibble::as_tibble(rownames = "GeneID") %>%
# dplyr::mutate(contrast = "48H vs 24H"),
# Fs_ddsObj.wald_DEGlist[[2]] %>%
# tibble::as_tibble(rownames = "GeneID") %>%
# dplyr::mutate(contrast = "1w vs 48H"),
# Fs_ddsObj.wald_DEGlist[[3]] %>%
# tibble::as_tibble(rownames = "GeneID") %>%
# dplyr::mutate(contrast = "3w vs 1w"),
# Fs_ddsObj.wald_DEGlist[[4]] %>%
# tibble::as_tibble(rownames = "GeneID") %>%
# dplyr::mutate(contrast = "4w vs 3w"),
# )
Get session info.
devtools::session_info()
## ─ Session info ───────────────────────────────────────────────────────────────
## setting value
## version R version 4.2.2 (2022-10-31)
## os macOS Big Sur ... 10.16
## system x86_64, darwin17.0
## ui X11
## language (EN)
## collate en_US.UTF-8
## ctype en_US.UTF-8
## tz Europe/Berlin
## date 2023-11-22
## pandoc 3.1.6.2 @ /usr/local/bin/ (via rmarkdown)
##
## ─ Packages ───────────────────────────────────────────────────────────────────
## package * version date (UTC) lib source
## annotate 1.76.0 2022-11-01 [1] Bioconductor
## AnnotationDbi 1.60.2 2023-03-10 [1] Bioconductor
## Biobase * 2.58.0 2022-11-01 [1] Bioconductor
## BiocGenerics * 0.44.0 2022-11-01 [1] Bioconductor
## BiocParallel 1.32.6 2023-03-17 [1] Bioconductor
## Biostrings 2.66.0 2022-11-01 [1] Bioconductor
## bit 4.0.5 2022-11-15 [1] CRAN (R 4.2.0)
## bit64 4.0.5 2020-08-30 [1] CRAN (R 4.2.0)
## bitops 1.0-7 2021-04-24 [1] CRAN (R 4.2.0)
## blob 1.2.4 2023-03-17 [1] CRAN (R 4.2.0)
## bslib 0.5.1 2023-08-11 [1] CRAN (R 4.2.0)
## cachem 1.0.8 2023-05-01 [1] CRAN (R 4.2.0)
## callr 3.7.3 2022-11-02 [1] CRAN (R 4.2.0)
## cli 3.6.1 2023-03-23 [1] CRAN (R 4.2.0)
## codetools 0.2-19 2023-02-01 [1] CRAN (R 4.2.0)
## colorspace 2.1-0 2023-01-23 [1] CRAN (R 4.2.0)
## crayon 1.5.2 2022-09-29 [1] CRAN (R 4.2.0)
## DBI 1.1.3 2022-06-18 [1] CRAN (R 4.2.0)
## DelayedArray 0.24.0 2022-11-01 [1] Bioconductor
## DESeq2 * 1.38.3 2023-01-19 [1] Bioconductor
## devtools 2.4.5 2022-10-11 [1] CRAN (R 4.2.0)
## digest 0.6.33 2023-07-07 [1] CRAN (R 4.2.0)
## dplyr * 1.1.3 2023-09-03 [1] CRAN (R 4.2.0)
## ellipsis 0.3.2 2021-04-29 [1] CRAN (R 4.2.0)
## evaluate 0.22 2023-09-29 [1] CRAN (R 4.2.2)
## fansi 1.0.5 2023-10-08 [1] CRAN (R 4.2.2)
## farver 2.1.1 2022-07-06 [1] CRAN (R 4.2.0)
## fastmap 1.1.1 2023-02-24 [1] CRAN (R 4.2.0)
## fdrtool 1.2.17 2021-11-13 [1] CRAN (R 4.2.0)
## forcats * 1.0.0 2023-01-29 [1] CRAN (R 4.2.0)
## foreach 1.5.2 2022-02-02 [1] CRAN (R 4.2.0)
## fs 1.6.3 2023-07-20 [1] CRAN (R 4.2.0)
## geneplotter 1.76.0 2022-11-01 [1] Bioconductor
## generics 0.1.3 2022-07-05 [1] CRAN (R 4.2.0)
## GenomeInfoDb * 1.34.9 2023-02-02 [1] Bioconductor
## GenomeInfoDbData 1.2.9 2023-05-04 [1] Bioconductor
## GenomicRanges * 1.50.2 2022-12-16 [1] Bioconductor
## ggplot2 * 3.4.4 2023-10-12 [1] CRAN (R 4.2.2)
## ggsci 3.0.0 2023-03-08 [1] CRAN (R 4.2.0)
## glue 1.6.2 2022-02-24 [1] CRAN (R 4.2.0)
## gtable 0.3.4 2023-08-21 [1] CRAN (R 4.2.0)
## hms 1.1.3 2023-03-21 [1] CRAN (R 4.2.0)
## htmltools 0.5.6.1 2023-10-06 [1] CRAN (R 4.2.2)
## htmlwidgets 1.6.2 2023-03-17 [1] CRAN (R 4.2.0)
## httpuv 1.6.11 2023-05-11 [1] CRAN (R 4.2.2)
## httr 1.4.7 2023-08-15 [1] CRAN (R 4.2.0)
## IHW * 1.26.0 2022-11-01 [1] Bioconductor
## IRanges * 2.32.0 2022-11-01 [1] Bioconductor
## iterators 1.0.14 2022-02-05 [1] CRAN (R 4.2.0)
## jquerylib 0.1.4 2021-04-26 [1] CRAN (R 4.2.0)
## jsonlite 1.8.7 2023-06-29 [1] CRAN (R 4.2.0)
## KEGGREST 1.38.0 2022-11-01 [1] Bioconductor
## knitr 1.44 2023-09-11 [1] CRAN (R 4.2.0)
## labeling 0.4.3 2023-08-29 [1] CRAN (R 4.2.0)
## later 1.3.1 2023-05-02 [1] CRAN (R 4.2.2)
## lattice 0.21-9 2023-10-01 [1] CRAN (R 4.2.2)
## lifecycle 1.0.3 2022-10-07 [1] CRAN (R 4.2.0)
## locfit 1.5-9.8 2023-06-11 [1] CRAN (R 4.2.0)
## lpsymphony 1.26.3 2023-01-19 [1] Bioconductor (R 4.2.2)
## lubridate * 1.9.3 2023-09-27 [1] CRAN (R 4.2.0)
## magrittr 2.0.3 2022-03-30 [1] CRAN (R 4.2.0)
## Matrix 1.5-4.1 2023-05-18 [1] CRAN (R 4.2.0)
## MatrixGenerics * 1.10.0 2022-11-01 [1] Bioconductor
## matrixStats * 1.0.0 2023-06-02 [1] CRAN (R 4.2.0)
## memoise 2.0.1 2021-11-26 [1] CRAN (R 4.2.0)
## mime 0.12 2021-09-28 [1] CRAN (R 4.2.0)
## miniUI 0.1.1.1 2018-05-18 [1] CRAN (R 4.2.0)
## munsell 0.5.0 2018-06-12 [1] CRAN (R 4.2.0)
## myTAI * 1.0.1.9000 2023-10-03 [1] Github (drostlab/myTAI@e159136)
## pillar 1.9.0 2023-03-22 [1] CRAN (R 4.2.0)
## pkgbuild 1.4.2 2023-06-26 [1] CRAN (R 4.2.0)
## pkgconfig 2.0.3 2019-09-22 [1] CRAN (R 4.2.0)
## pkgload 1.3.3 2023-09-22 [1] CRAN (R 4.2.0)
## png 0.1-8 2022-11-29 [1] CRAN (R 4.2.0)
## prettyunits 1.2.0 2023-09-24 [1] CRAN (R 4.2.0)
## processx 3.8.2 2023-06-30 [1] CRAN (R 4.2.0)
## profvis 0.3.8 2023-05-02 [1] CRAN (R 4.2.0)
## promises 1.2.1 2023-08-10 [1] CRAN (R 4.2.2)
## ps 1.7.5 2023-04-18 [1] CRAN (R 4.2.0)
## purrr * 1.0.2 2023-08-10 [1] CRAN (R 4.2.2)
## R6 2.5.1 2021-08-19 [1] CRAN (R 4.2.0)
## RColorBrewer 1.1-3 2022-04-03 [1] CRAN (R 4.2.0)
## Rcpp 1.0.11 2023-07-06 [1] CRAN (R 4.2.0)
## RCurl 1.98-1.12 2023-03-27 [1] CRAN (R 4.2.0)
## readr * 2.1.4 2023-02-10 [1] CRAN (R 4.2.0)
## remotes 2.4.2.1 2023-07-18 [1] CRAN (R 4.2.2)
## rlang 1.1.1 2023-04-28 [1] CRAN (R 4.2.0)
## rmarkdown 2.25 2023-09-18 [1] CRAN (R 4.2.2)
## RSQLite 2.3.1 2023-04-03 [1] CRAN (R 4.2.0)
## rstudioapi 0.15.0 2023-07-07 [1] CRAN (R 4.2.0)
## S4Vectors * 0.36.2 2023-02-26 [1] Bioconductor
## sass 0.4.7 2023-07-15 [1] CRAN (R 4.2.0)
## scales 1.2.1 2022-08-20 [1] CRAN (R 4.2.0)
## sessioninfo 1.2.2 2021-12-06 [1] CRAN (R 4.2.0)
## shiny 1.7.5.1 2023-10-14 [1] CRAN (R 4.2.2)
## slam 0.1-50 2022-01-08 [1] CRAN (R 4.2.0)
## stringi 1.7.12 2023-01-11 [1] CRAN (R 4.2.0)
## stringr * 1.5.0 2022-12-02 [1] CRAN (R 4.2.0)
## SummarizedExperiment * 1.28.0 2022-11-01 [1] Bioconductor
## tibble * 3.2.1 2023-03-20 [1] CRAN (R 4.2.0)
## tidyr * 1.3.0 2023-01-24 [1] CRAN (R 4.2.0)
## tidyselect 1.2.0 2022-10-10 [1] CRAN (R 4.2.0)
## tidyverse * 2.0.0 2023-02-22 [1] CRAN (R 4.2.0)
## timechange 0.2.0 2023-01-11 [1] CRAN (R 4.2.0)
## tzdb 0.4.0 2023-05-12 [1] CRAN (R 4.2.2)
## urlchecker 1.0.1 2021-11-30 [1] CRAN (R 4.2.0)
## usethis 2.2.2 2023-07-06 [1] CRAN (R 4.2.0)
## utf8 1.2.3 2023-01-31 [1] CRAN (R 4.2.0)
## vctrs 0.6.4 2023-10-12 [1] CRAN (R 4.2.2)
## vroom 1.6.4 2023-10-02 [1] CRAN (R 4.2.2)
## withr 2.5.1 2023-09-26 [1] CRAN (R 4.2.0)
## xfun 0.40 2023-08-09 [1] CRAN (R 4.2.2)
## XML 3.99-0.14 2023-03-19 [1] CRAN (R 4.2.0)
## xtable 1.8-4 2019-04-21 [1] CRAN (R 4.2.0)
## XVector 0.38.0 2022-11-01 [1] Bioconductor
## yaml 2.3.7 2023-01-23 [1] CRAN (R 4.2.0)
## zlibbioc 1.44.0 2022-11-01 [1] Bioconductor
##
## [1] /Library/Frameworks/R.framework/Versions/4.2/Resources/library
##
## ──────────────────────────────────────────────────────────────────────────────